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Abstract 

A dispersion-managed optical system with step-wise periodical variation of dis- 
persion is studied in a strong dispersion map limit in the framework of path-averaged 
Gabitov-Turitsyn equation. The soliton solution is obtained by iterating the path- 
averaged equation analytically and numerically. An efficient numerical algorithm 
for obtaining of DM soliton shape is developed. The envelope of soliton oscillating 
tails is found to decay exponentially in time while the oscillations are described by a 
quadratic law. 
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A dispersion-managedS (DM) optical fiber is designed to create a low (or even zero) 
path-averaged dispersion by periodically alternating dispersion sign along an optical line 
which dramatically reduces pulse broadening. Recently dispersion management has become 
an essential technology for development of ultrafast high-bit-rate optical communication 
linesH. Lossless propagation of optical pulse in DM fiber is described by a nonlinear 
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Schrodinger equation (NLS) with periodically varying dispersion d(z): 

iu z + d(z)utt + \u\ 2 u = 0, (1) 

where u is the envelope of optical pulse, z is the propagation distance and all quantities 
are made dimensionless. Consider a two-step periodic dispersion map: d(z) = (d) + d(z), 
where d(z) = d% for < z + nL < Li and d(z) = d 2 for L\ < z + nL < L, L = Li + L 2 
is a dispersion map period, (d) is the path-averaged dispersion, di,d 2 are the amplitudes of 
dispersion variation subjected to condition d\L\ + d 2 L 2 = and n is an arbitrary integer 
number. 

A nonlinearity can be treated as a small perturbation on a scales of dispersion map 
period L provided a characteristic nonlinear length Z n \ of the pulse is large: Z n \ 3> L, 
where Z n i = l/\p\ 2 and p is a typical pulse amplitude. Then the eq. flU) is reduced to a 
path-averaged Gabitov-TuritsynS model: 

itp z (u) -cu 2 (d)^ + R(^,to) = 0, (2) 
where 

• sA 

sin 



xi/;*(uj 3 )S(uji + uj 2 — cl> 3 — u)cLuicLu 2 duj3, (3) 

a 9999 rr-i- • i?« * tj2 f d(z')dz' . 

A = + uo 2 — uo^ — Lu 1 , s = ctiLi is a dispersion map strength, ip = ue L i /2 is a 

slow function of z on a scale L and V>(u;) — Xf^o ipifye^dt is a Fourier component of ip. The 
Gabitov-Turitsyn model is well supported by numerical simulations!'!. 

Consider DM soliton solution ip = A(t)e lXz (A is real) of the Gabitov-Turitsyn Eq. (Q) 
then returning to t-space one gets: 

- \A + (d)A tt = ^-f Ci(p^)A(t x +t)x 



2ns J s 

Afo + ^Aih + h + t^dh, (4) 
where Ci(x) = J^cosx/xdx. It was found numerically!! that the Gaussian ansatz 



A GauS s = pexp ( - ^t 2 ), (5) 

where p, (3 are real constants, is a rather good approximation for the DM soliton solution. 
Thus the Eq. (|5|) can be effectively used as zero approximation for solving Eq. (|J) by itera- 
tions which was done in Ref. for (d) = 0. Following^ one can easily make a generalization 
for the case of small but nonzero average dispersion \d \ <C \di\ and obtain a set of two 
transcendental Eqs. from the series expansion of a first iteration in powers of t 2 : 

A = —f3{d) + (arcsinh——^ h c.c), 

A = -313(d) 



2p 2 ( I s + t V3 . , 3s - 1 \ 
H W 1 arcsinh h c.c. , 

where c.c. means complex conjugation. The Eqs. @ determine the parameters /3, p of the 
Gaussian ansatz @ as a functions of system parameters A, s. 

The solution of Eq. @ was also obtained by means of iterating this Eq. numerically. 
n + 1th iteration A^ n+1 ' is given by: 

where the functional R(A,uj) is defined in (|3|), Q n is a stabilizing factor given by 



r \\+\{d)\^) 

and F~ x is a Backward Fourier transform. This numerical iteration scheme was also used in 
Ref.0 except that a Petviashvili stabilizing factoid was used there instead of Q n . But both 
stabilizing factors results in the convergence of iteration scheme to the same solution of (||). 

The main obstacle in numerical iteration scheme (||), (|7|) is the computation of integral 
term R(A, uj) which generally require iV 3 operations for each iteration, where N is a number 
of grid points in uj or t-space. Here we introduce much more efficient numerical algorithm 
for calculation of R(A,u). 

Rewriting the kernel of R(ip,uj) via parametric integral: 



and using definition of A one gets from (|3]) : 

R(A,u) = — ^- f S/2 ds'e~^ 2 f i^M 

S(Z7T) Z J-s/2 J 

x A {sl) (uj 2 )A {sI) * {u)$8{u)\ + lu 2 -u 3 - u)du 1 du 2 du 3 , (9) 

where A^ s '\uj) = A(uj)e ls,uj2 . In t-space this expression takes the form 

F-\R{A,uj)) = - f 2 ds'G {s,) {¥ s '\t)), (10) 

S J-s/2 

where ^^(t) = |A( s ')(t)| 2 A( s ')(t) and G( s ') is an integral operator corresponding to a multi- 
plication operator G^^^'^o;)) = e ls,uj2 ^( s '} in cu-space. It follows from the Eqs. (Q), (p!U[) 
that numerical procedure for calculation of R(A, u) includes four steps: 

(i) The Backward Fourier Transform of A^(u) = A(u))e %s>u)2 for every value of s' . 

(ii) A calculation of ^( s ')(t) from A^ s '\t). 

(hi) The Forward Fourier Transform of ^ s '(t). 

(iv) A numerical integration (summation) of e^'^^f^'^u) over s' for every value of cu. 

The Forward and Backward Fourier Transforms were performed with the Fast Fourier 
Transform which requires NLog 2 (N) operations. A total number of operations for one 
iteration is about 4MNLog 2 (N), where M is a number of grid points for integration over 
s'. We used the following typical values for numerical solution of @: N = 8192; M = 
800. One iteration on Alpha 500MHz workstation requires about 30 seconds for 16-bytes 
(32 digits) precision. Thus numerical scheme (i)-(iv) provides dramatic improvement of 
numerical performance. 8192 3 operations would takes 30 days on the same workstation. 
Note that the proposed efficient numerical algorithm can be generalized to include optical 
fiber losses and amplifiers. 

Fig.la,b show the dependence of a root mean square pulse width Trms = 
yj J t 2 A 2 dt/ J A 2 dt on a quasimomentum A obtained from (i) the first iteration of the Eq. (|]) 
using values of (3,p resulting from the Eqs. (|6|) (dotted curves); (ii) a variational approach 



(see e.g. Eqs. (13), (14) in Ref.0) represented by dashed lines; (iii) a full numerical solu- 
tion of the Eq. (f|) (solid lines). The explicit expression T RMS = l/y/2/3 for the Gaussian 
pulse shape is used for calculation of dashed curves. The solid curve is shown in Fig. lb 
only for upper branch I of solution because numerical iteration scheme for negative average 
dispersion (d) = —0.01 diverges on lower branch II which is in agreement with Ref.0. We 
also calculated a time-averaged optical power P = J A 2 dt and found that P(X) dependence 
following from the first iteration and the variational approach Ref.0 reproduce a full nu- 
merical solution of the Eq. @ with high accuracy (~ 1%). One can conclude that both the 
Eqs. (^|) and the variational approach.0 predict P(X) with a high accuracy while Trms(X) 
dependence is reproduced by the first iteration of the Eq. @ with better accuracy (~ 2%) 
compare with accuracy (~ 40%) of the variational approach. 

There is an essential difference of our numerical simulation in comparison with numerical 
results of Ref.0 concerning upper branch I for the negative average dispersion. After about 
50 iterations of the Eq. (^) a numerical instability was detected on the tails of DM soliton for 
(d) = —0.01. (Presumably this instability was not found in Ref.0 because a few iterations 
was considered there). A finer numerical grid slows down numerical instability growth but 
does not kill it. The instability slows down as (d) — > and for (d) > there is no numerical 
instability. Thus the solid curve in Fig. lb for (d) = —0.01 can only be formally attributed 
to DM soliton and the question about existence of DM soliton for the negative average 
dispersion is still open. It is possible that the instability within the numerical iteration 
scheme results from a resonance of DM soliton tails with linear waves0. But there is another 
alternative that DM soliton solution does not exist for any negative average dispersion value 
and instead of DM soliton on can observe a long-lived quasi-stable structure. Note that the 
existence of DM soliton for nonnegative average dispersion for the Eq. (^) was proved in 
Ref.fi In addition it was proved in Ref." that even if DM soliton exists for (d) < it can 
not realize a minimum of the Hamiltonian of the Eq. (|2[) for fixed P. This can indicate 
that DM soliton is unstable in that case. Related result@ is the nonexistence criterion for a 
periodic solution of the Eq. (fj) for a negative enough average dispersion. But Refs.EUli do 
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not give any statement about the existence of DM soliton for small negative (d) . 

Fig. 2 shows a typical shape of DM soliton. This is the first to the best of my knowledge 
high precision numerical solution of the Eq. Note that solid curve dips do not reach 
t— axes only because of finite size of numerical grid. On can conjecture from Fig. 2 the 
asymptotic of DM soliton is given by 

A asymp (t) = f(t) cos (t 2 [a + a(t)}) exp (-b\t\), (11) 

where ao, b are constants and f(t)/\t\, a(t) are slow functions oft. An analysis of fast oscilla- 
tions in integral term of the Eq. (Q) allows to show that f(t) = c\t\+0(l), ao = l/2s, a(t) = 
ai/|t|+a 2 A 2 + 0(l/|t| 3 ) for \t\ -> oo, (d) -> 0, where c, ai, a,2 are constants. Dashed curve in 
Fig. 2 shows A 2 asymp (t) dependence for c = 11.9654, b = 3.04515, a x = 1.41364, a 2 = 1.51023 
which is in very good agreement with asymptotic of numerical solution of the Eq. (solid 
curve). Thus the envelope of DM soliton oscillating tails decays exponentially while the os- 
cillations are described by a quadratic law. Detailed consideration of the asymptotic solution 
is outside the scope of this Letter. 

The author thanks M. Chertkov, I.R. Gabitov, E.A. Kuznetsov and V. Zharnitsky for 
helpful discussions. 
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Figure captions: 



Fig.l. T RMS for s = l,(d) = 0.01 (a) and (d) = -0.01 (b). Branches I, II for (d) = -0.01 
correspond to two branches of analytical solution. 

Fig.2. DM soliton shape (curve 1) versus Eq. (11) (curve 2) for (d) — 0, s — 1, A = 1. A(t) 
is an even function. 
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Fig. 1a 
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Fig. 1b 
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